package PlotPourri;

import Jama.Matrix;
import java.io.File;
import java.io.FileNotFoundException;
import java.util.ArrayList;
import java.util.List;
import java.util.Scanner;

public class Observation extends ArrayList<double[]> {

    private String path = "";
    private double period, distance, reddening = 0;
    private int name, magnitudePosition, errorPosition, julianDatePosition = 0;
    private int phasePosition = -1;
    private BandName band = BandName.values()[0];

    public Observation() {
    }

    public Observation(String path, int name) {
        this.path = path;
        this.name = name;
        parse(new File(path));
    }

    public Observation(File file) {
        this.path = file.getPath();
        parse(file);
    }

    public double getPhase(int index) {
        if (getPhasePosition() > 0) return get(index)[getPhasePosition()];
        return ((getJulianDate(index) / getPeriod()) - (int) (getJulianDate(index) / getPeriod()));
    }

    public double getMagnitude(int index) {
        return get(index)[getMagnitudePosition()];
    }

    public double getError(int index) {
        return get(index)[getErrorPosition()];
    }

    public double getJulianDate(int index) {
        return get(index)[getJulianDatePosition()];
    }

    private void parse(File file) {
        Scanner sc = null;
        try {
            sc = new Scanner(file);
        } catch (NullPointerException ex) {
            System.err.print(ex + "\n");
        } catch (FileNotFoundException ex) {
            System.err.print(ex + "\n");
        }
        if (sc != null) {
            while (sc.hasNextLine()) {
                String line = sc.nextLine();
                while (line.startsWith(" ")) line = line.substring(1);
                double[] dbls = new double[line.split("\\s+").length];
                for (int i = 0; i < line.split("\\s+").length; i++) dbls[i] = Double.parseDouble(line.split("\\s+")[i]);
                add(dbls);
            }
            sc.close();
        }
    }

    public double[][] getData(int order, int xsize) {
        return getData(getCoefficients(order), order, xsize);
    }

    public double[][] getData(double[] coefficients, int order, int xsize) {
        double Ao = getAo(coefficients);
        List<Double> A = getA(coefficients);
        List<Double> phi = getPhi(coefficients);
        double[][] data = new double[xsize + 1][2];
        for (int i = 0; i <= xsize; i++) {
            double phase = (double) i / xsize;
            double sum = 0;
            for (int k = 0; k < order; k++) sum += (double) A.get(k) * Math.sin((k + 1) * 2 * Math.PI * phase + phi.get(k));
            double v = Ao + sum;
            data[i][0] = phase;
            data[i][1] = v;
        }
        return data;
    }

    public double maxLight(int order, int xsize) {
        double[][] data = getData(order, xsize);
        double smallest = Double.POSITIVE_INFINITY;
        for (double[] array : data) if (array[1] < smallest) smallest = array[1];
        return smallest;
    }

    public double maxLight(double[] coefficients, int order, int xsize) {
        double[][] data = getData(coefficients, order, xsize);
        double smallest = Double.POSITIVE_INFINITY;
        double position = 0;
        for (double[] array : data) if (array[1] < smallest) {
                smallest = array[1];
                position = array[0];
            }
        return position;
    }

    public double[][] getShiftedData(double[] coefficients, int order, int xsize) {
        double Ao = getAo(coefficients);
        List<Double> A = getA(coefficients);
        List<Double> phi = getPhi(coefficients);
        double[][] data = getData(coefficients, order, xsize);
        double maximum = Double.POSITIVE_INFINITY;
        double Po = 0;
        for (int i = 0; i < data.length; i++) if (data[i][1] < maximum) {
                maximum = data[i][1];
                Po = data[i][0];
            }
        for (int i = 0; i <= xsize; i++) {
            double phase = (double) i / xsize;
            double sum = 0;
            for (int k = 0; k < order; k++)
                sum += (double) A.get(k) * Math.sin((k + 1) * 2 * Math.PI * phase + phi.get(k) + (2 * Math.PI * Po * (k + 1)));
            double v = Ao + sum;
            data[i][0] = phase;// + Po % 100;
            data[i][1] = v;
        }
        return data;
    }

    public double[][] getShiftedData(int order, int xsize) {
        return getShiftedData(getCoefficients(order), order, xsize);
    }

    public double[] getAkBk(int order, double threshold) {
        double[] coefficients = getCoefficients(order);
        List<Double> A = getA(coefficients);
        for (Double d : A) if (Math.abs(d) > threshold) {
                double[] zeros = new double[order * 2];
                for (int i = 0; i < zeros.length; i++) zeros[i] = 0;
                return zeros;
            }
        List<Double> phi = getPhi(coefficients);
        double[] a = new double[A.size()];
        double[] b = new double[A.size()];
        double phmax = maxLight(coefficients, order, 100);
        for (int k = 0; k < A.size(); k++) {
            a[k] = A.get(k) * Math.cos(phi.get(k) + 2 * Math.PI * k * phmax) - coefficients[0]; // a_k
            b[k] = -A.get(k) * Math.sin(phi.get(k) + 2 * Math.PI * k * phmax) - coefficients[0]; // b_k
        }
        double[] result = new double[A.size() * 2];
        for (int i = 0; i < a.length; i++) {
            result[i * 2] = a[i];
            result[i * 2 + 1] = b[i];
        }
        return result;
    }

    public double getAo(int order) {
        return getAo(getCoefficients(order));
    }

    private double getAo(double[] coefficients) {
        return coefficients[0];// - getReddening();
    }

    public double getV(int order, double phase) {
        double[] coefficients = getCoefficients(order);
        double Ao = coefficients[0];
        List<Double> A = getA(coefficients);
        List<Double> phi = getPhi(coefficients);
        return getV(order, phase, Ao, A, phi);
    }

    private double getV(int order, double phase, double Ao, List<Double> A, List<Double> phi) {
        double sigma = 0;
        for (int k = 0; k < order; k++) sigma += (double) A.get(k) * Math.sin((k + 1) * 2 * Math.PI * phase + phi.get(k));
        return Ao + sigma;
    }

    List<Double> getA(double[] coefficients) {
        List<Double> A = new ArrayList<Double>((coefficients.length - 1) / 2);
        for (int i = 1; i < coefficients.length && i < (coefficients.length - 1) / 2 + 1; i++) A.add(coefficients[i]);
        return A;
    }

    List<Double> getPhi(double[] coefficients) {
        List<Double> Phi = new ArrayList<Double>((coefficients.length - 1) / 2);
        for (int i = (coefficients.length - 1) / 2 + 1; i < coefficients.length; i++) Phi.add(coefficients[i]);
        return Phi;
    }

    public Observation quickSort() {
        Observation less = new Observation();
        Observation greater = new Observation();
        if (size() <= 1) return this;
        double[] pivot = get(size() / 2);
        for (double[] x : this)
            if (x == pivot) continue;
            else if (x[0] <= pivot[0]) less.add(x);
            else greater.add(x);
        Observation concatenated = new Observation();
        concatenated.setName(getName());
        concatenated.setPeriod(getPeriod());
        for (double[] x : less.quickSort()) concatenated.add(x);
        concatenated.add(pivot);
        for (double[] x : greater.quickSort()) concatenated.add(x);
        return concatenated;
    }

    public double getDisp(int order) {
        Observation sortedData = quickSort();
        double[] coefficients = sortedData.getCoefficients(order);
        double Ao = coefficients[0];
        List<Double> A = getA(coefficients);
        List<Double> phi = getPhi(coefficients);
        List<Double> vValues = new ArrayList<Double>();
        for (int i = 0; i < sortedData.size(); i++)
            vValues.add(Math.pow((sortedData.getMagnitude(i)) - getV(order, sortedData.getPhase(i), Ao, A, phi), 2));
        double disp = 0;
        for (double d : vValues) disp += d;
        return disp / vValues.size();
    }

    public double getDisp(int order, int xsize, double threshold) {
        Observation sortedData = quickSort();
        double[] coefficients = sortedData.getCoefficients(order);
        if (hasOscillations(sortedData, order, xsize)) return Double.POSITIVE_INFINITY;
        double Ao = coefficients[0];
        List<Double> A = getA(coefficients);
        List<Double> phi = getPhi(coefficients);
        List<Double> vValues = new ArrayList<Double>();
        for (int i = 0; i < sortedData.size(); i++)
            vValues.add(Math.pow((sortedData.getMagnitude(i)) - getV(order, sortedData.getPhase(i), Ao, A, phi), 2));
        double disp = 0;
        for (double d : vValues) disp += d;
        System.out.println("Dispersion at order " + order + ": " + (disp / vValues.size()));
        return disp / vValues.size();
    }

    public boolean hasOscillations(Observation sortedData, int order, int xsize) {
        double[][] data = getData(order, xsize);
        for (int i = 0; i < sortedData.size() - 1; i++) {
            if (i < sortedData.getPhase(i)) continue;
            if (i > sortedData.getPhase(i + 1)) break;
            double previous = 0;
            boolean increasing = false;
            boolean decreasing = false;
            for (double[] d : data) {
                if (d[1] > previous) increasing = true;
                if (d[1] < previous) decreasing = true;
                if (increasing && decreasing) return true;
                previous = d[1];
            }
        }
        return false;
    }

    public double[] getCoefficients(int order) {
        int fit = 2 * order + 1;
        int dataSize = size();
        double[][] helper = new double[fit][dataSize];
        for (int i = 0; i < dataSize; i++) helper[0][i] = 1;
        for (int i = 1; i < fit; i++)
            for (int j = 0; j < dataSize; j++) {
                double t = 2 * Math.PI * getPhase(j);
                if (i <= order) helper[i][j] = Math.sin(i * t);
                else helper[i][j] = Math.cos((i - order) * t);
            }
        double[] B = new double[fit];
        for (int i = 0; i < fit; i++) {
            double sum = 0;
            for (int j = 0; j < dataSize; j++) sum += getMagnitude(j) * helper[i][j];
            B[i] = sum;
        }
        // get our matrix
        Matrix newHelper = new Matrix(helper);
        Matrix newA = newHelper.times(newHelper.transpose());
        Matrix matrixB = new Matrix(B, fit);
        Matrix newX;
        try {
            newX = newA.solve(matrixB);
        } catch (RuntimeException ex) {
            double[] c = new double[fit];
            for (int i = 0; i < fit; i++) c[i] = 0;
            return c;
        }
        // take our A0, ai's, bi's and change to A's phi's
        double[] ourAns = new double[fit];
        ourAns[0] = newX.get(0, 0);
        for (int i = 1; i < order + 1; i++) ourAns[i] = Math.atan((newX.get(i + order, 0)) / (newX.get(i, 0)));
        for (int i = (order + 1); i < fit; i++) ourAns[i] = newX.get(i - order, 0) / (Math.cos(ourAns[i - order]));
        //fix the order, A's before phi's
        double[] c = new double[fit];
        c[0] = ourAns[0];
        for (int i = 1; i < (order + 1); i++) c[i] = ourAns[i + order];
        for (int i = (order + 1); i < fit; i++) c[i] = ourAns[i - order];
        //fix it so A's are positive
        for (int i = 1; i < (order + 1); i++)
            if (c[i] < 0) {
                c[i + order] = c[i + order] + Math.PI;
                c[i] = -1.0 * c[i];
            }
        //fix it so phis are positive
        for (int i = (order + 1); i < fit; i++) if (c[i] < 0) c[i] = c[i] + 2 * Math.PI;
        return c;
    }

    public double getPeriod() {
        return period;
    }

    public void setPeriod(double period) {
        this.period = period;
    }

    public double getDistance() {
        return distance;
    }

    public void setDistance(double distance) {
        this.distance = distance;
    }

    public double getReddening() {
        return reddening;
    }

    public void setReddening(double reddening) {
        this.reddening = reddening;
    }

    public String getPath() {
        return path;
    }

    public void setPath(String path) {
        this.path = path;
    }

    public int getName() {
        return name;
    }

    public void setName(int name) {
        this.name = name;
    }

    public BandName getBand() {
        return band;
    }

    public void setBand(BandName band) {
        this.band = band;
    }

    public int getPhasePosition() {
        return phasePosition;
    }

    public void setPhasePosition(int phasePosition) {
        this.phasePosition = phasePosition;
    }

    public int getMagnitudePosition() {
        return magnitudePosition;
    }

    public void setMagnitudePosition(int magnitudePosition) {
        this.magnitudePosition = magnitudePosition;
    }

    public int getErrorPosition() {
        return errorPosition;
    }

    public void setErrorPosition(int errorPosition) {
        this.errorPosition = errorPosition;
    }

    public int getJulianDatePosition() {
        return julianDatePosition;
    }

    public void setJulianDatePosition(int julianDatePosition) {
        this.julianDatePosition = julianDatePosition;
    }
}
